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In this paper we discuss some general aspects of the gravitational wave background arising from 
post-inflationary short-lasting cosmological events such as phase transitions. We concentrate on the 
physics which determines the shape and the peak frequency of the gravitational wave spectrum. We 
then apply our general findings to the case of bubble collisions during a first order phase transition 
and compare different results in the recent literature. 
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I. INTRODUCTION 

In cosmology there are several situations in which a 
stochastic gravitational wave (GW) background can be 
generated. For example, inflation leads to the quantum 
generation of gravitons which are relevant at very large 
wavelength. Here we are interested in gravitational waves 
produced after inflation, e.g. during preheating ll[ or 
during the electroweak phase transition [2|, y, |4|, [5|, la] • In 
these situations the gravitational waves are sourced by a 
transverse (tensor type) anisotropic stress in the cosmic 
fluid. As these stresses are generated causally after in- 
flation, they have a finite correlation length R which is 
limited by the Hubble scale. In the cases we want to dis- 
cuss in this work, the anisotropic stress is non-vanishing 
for a finite duration fi^ 1 which is assumed to be smaller 
than the Hubble time, (3 ^S> H.(t*). Here denotes the 
(conformal) time when the phase transition (or preheat- 
ing) begins and it ends at + The time scale 
and the correlation length are related by some velocity 
v < 1, R ~ v/{3. In the literature, the peak of the energy 
spectrum of the GWs has been found both at wavenum- 
ber k ~ (3 0, i, H or k ~ i?" 1 @, 0], and the question 
of the correct peak frequency of the GW spectrum from 
cosmological sources is still under debate [8( . While it is 
not contested that causality implies that the GW spec- 
trum scales as "fiogl oc ^ 3 ^ or sma H frequencies, k < (3, 
it is still unclear what precisely determines the position 
of the peak and how the GW power spectrum decays for 
large frequencies. 

In this paper we want to address these questions. We 
shall clearly identify the properties of the anisotropic 
stress which determine the peak position and the decay 
law at large frequencies and we shall clarify several spe- 
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cific examples. 

In the next section we relate the gravity wave energy 
spectrum to the diagonal of the anisotropic stress spec- 
trum. In Section III we discuss several possibilities for 
the unequal time correlator of the stress tensor and de- 
termine the resulting peak frequency. In Section IV we 
study in detail the case of bubble collisions which has 
been discussed in two recent papers [H, Q with conflicting 
results. We clarify the difference of the two treatments 
and argue that an unphysical assumption in [|| (discon- 
tinuity of the anisotropic stress at the end of the transi- 
tion) leads to a peak position that is not at fc pca k ~ /3, 
as found in We also reveal the origin of the mild 1/fc 
decay of the spectrum obtained in [5| and show that it is 
quite fragile to small modifications in the modeling. In 
Section V we conclude. 

Notation: we work in conformal time called t, so that 
the perturbed metric is given by 

ds 2 = a 2 (t) (-dt 2 + (6ij + h^dx 1 dx j ) 

where hij is transverse traceless, i.e. a gravitational wave 
perturbation. We define the conformal Hubble parameter 
TC=^/a=a/a = aH. The scale factor is normalized 
to unity today, so that conformal wavenumber k becomes 
the physical wavenumber/frequency today. 



II. GW ENERGY DENSITY SPECTRUM FROM 
A STOCHASTIC SHORT-LASTING SOURCE 

We consider a gravitational wave source, i.e., a tensor 
type (transverse) anisotropic stress coming from either 
colliding bubbles, or turbulence, or a stochastic scalar or 
vector {e.g. magnetic) field, etc., Hij(x,t). This leads to 
the generation of gravitational waves via the linearized 
Einstein equation 
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is the d'Alembert operator (in a cosmological back- 
ground), Hij is the dimensionless anisotropic stress and 
px the energy density of the source. We consider this 
source to be a statistically homogeneous and isotropic 
random variable with a power spectrum P s (k,t,t') de- 
fined by 

n ii (k,t)= 4( k ) n -4(M), (3) 

A=l,2 

(n A (k,t)U* B (k',t r )) = (2ir) 3 S(k-k>)6 AB P s (k,t,t') .(4) 

Here e^-(k) is a normalized polarization tensor (e.g. the 
hclicity basis) and we assume parity invariance so that 
both helicities have the same spectrum and are mutually 
uncorrelated. 

We want to consider short-lived sources. In the cos- 
mological context, a source is called short-lived if it is 
non-zero from some initial time £* until some final time 

+ 1/(3 with 1/(3 < H* 1 = li.- 1 ^*). In the short- 
lasting case, we can neglect the Hubble damping during 
the time when the source is active and we can write the 
wave equation in the form (we suppress the index A since 
the result is the same for both polarizations), a* = a(t*) 



(ff + mM) = 32 ^°- p * n(k,t) 



(5) 



At times t > t* + 1/(3 but still during the radiation era, 
the solution on sub-horizon scales, ft 3> Ti is 
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At times t < t* + l/(3 the integral in the above expres- 
sion only extends until t and the pre- factor a* /a can be 
neglected. Here 
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dt / dt'Psik^t'yt"*-"'*") . (9) 



The second term in ([7]), which is multiplied by e 2 , av- 
erages to zero over an oscillation period (see also [§}). 
Note that apart from the fact that the source is short 
lasting, we did not make any assumption about its time 
structure so far. 

The gravitational wave energy density is defined as 
Pgw{*) = (hij(x)h*j(x.))/(8irGa 2 ). Fourier transforming 
this expression and using h ~ k h we obtain 
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Here we assume that the gravitational wave is generated 
during the radiation dominated era. 

To determine the gravitational wave spectrum it suf- 
fices therefore to study 7i 2 k 3 Re[P s (k, k, ft)], which is a 
dimensionless quantity (note that ft 3 P s (k, t, t') is dimen- 
sionless, hence P s {k,uj,uj') has the dimension of time to 
the fifth power, we work in units with c = K = 1). 

Eq. (fTU|) is physically equivalent to Eq. 10.4.16 and 
following in Weinberg's book [9j, if one normalizes the 
latter to the critical energy density in the universe, ex- 
presses it per logarithmic unit of frequency, integrates it 
over directions and considers a stochastic source which 
is statistically homogeneous and isotropic such that 
(A^ iJm (k)i;*.(k,w)Tj m (k,w)) corresponds to P s (k,k,k). 

The result (fTU)) is very general for scales which enter 
the horizon during the radiation dominated era (i.e. fre- 
quencies larger than about 10 _11 Hz). We now analyze 
different physical situations and discuss the features of 
the expected gravitational wave spectrum. 



III. SOME GENERAL EXAMPLES 



is the time- Fourier transform of n(k, t). Gravitational 
waves are only sensitive to the diagonal of the Fourier 
transform of the anisotropic stress, \u)\ = k. 

The spectrum of the tensor perturbations at t > t* + 
1/(3, k -C H(t) becomes 

(h(k, t)h*(k',t)} 
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In this section we discuss four different forms for 
the unequal time power spectrum P s (k,t,t') defined in 
Eq. ([3]). These forms are quite general, and have been 
proposed already in Ref. [6J in the context of bubble col- 
lisions. We analyze the GW spectrum that arises in each 
of these cases, and we are mainly concerned with the 
time structure of P s (k,t,t'). To determine the gravita- 
tional wave spectrum we need the double time Fourier 
transform given in Eq. ([9]). We want to maintain sta- 
tistical homogeneity and isotropy of the source in space, 
because it is always justified in the cosmological context. 
From these properties it follows that the fc-dependence 
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of P s (k, k, k) due to the spatial structure of the source is 
given simply by the space Fourier transform of the source 
itself: statistical homogeneity and isotropy imply 

(n(x,t)n*(x',t')) = P s (|x-x'|,i,0 (ii) 

and therefore (with z = x — x', z = |z|) 

(n(k,f)n*(k',0) = 

= (2^) 3 £ 3 (k-k') f d S ze ikx P s (z,t,t'), (12) 

and 

P s (k,t,t') =4tt / dzz 2 S J^™t pJ z ,t,t'). (13) 
Jo fcz 

More specifically, for the illustrative purpose of this 
section, we assume that the anisotropic stress power spec- 
trum at equal time is separable, 



P s (k 1 t,t) = \F(k)\ 2 \g(t)f 



(14) 



If the source generating the gravitational waves is causal, 
this correlation function has compact support in space 
given by the correlation scale R, i.e. P s (\z\,t,t') = for 
|z| > R. Then its Fourier transform is analytic at k = 
which generically means that it is white noise on large 
scales (note that the tensor structure of the correlator 
can impose a different behavior, e.g. for magnetic fields, 
which have a k 2 spectrum on large scales [l(| )■ Fur- 
thermore, for the total energy in gravitational waves to 
remain finite, k 3 \F(k)\ 2 has to decay for k — + oo. A sim- 
ple Ansatz which satisfies these requirements and which 
has the correct dimensions is 



\F(k)f 



R :i 



i + {kpy 



(15) 



where R denotes the characteristic scale of the problem, 
typically the correlation scale. As we shall see in the ex- 
ample of colliding bubbles, it is more realistic to assume 
that R is time dependent, and this time dependence can 
affect the spectrum. The assumption of separability has 
some immediate consequences, namely that the slope of 
the spectrum changes at the frequency k ~ 1/R. With 
the above choice the change in slope is k . 

We relate the characteristic length scale R to the char- 
acteristic time scale (3 by a velocity v, R = v/(3. In the 
following, we analyze three forms for the function g(t). 
The first one is discontinuous 1 : 



9i(t) = 



1 U<t<t. 
else , 
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1 This time dependence is unphysical as it implies that the energy 
momentum tensor is discontinuous, but it is possible to have sit- 
uations where the energy momentum tensor changes very rapidly 
and which can therefore be approximated by a discontinuity. 



the second one is continuous but not differentiable at 
t = t* and t = t» + 1/(3, i.e. g(t) is in C° (but not in C 1 ) 
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and the third one is in C 1 (but not in C 2 ) at t» and 

i* + 1//3, 
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We now go on to analyze four different possibilities for 
the unequal time correlation function P s (k,t,t') which 
at equal times reproduce the form given in Eq. (|14[) with 
the functions F(k) and g(t) given above. 



A. Totally incoherent sources 

Let us first assume that the source at different times is 
not correlated, i.e. it is a sequence of very short events. 
We call such a source totally incoherent. In this case 

(n(k, t)n* (k, t')) = (2nfs(k - V) S Jtzfl Ps{kj t , t) . 

(19) 

We have introduced the time scale 1//3, the duration 
of the source, to take care of dimensions. For the 
anisotropic stress power spectrum we obtain 

Ps(k,t,t') = ^V(fc)| 2 | 5 (<)| 2 

P s (k,k,k) = r dt\g(t)\ 2 . 

In this situation, the spectrum P s (k, k, k) is not affected 
by the time Fourier transform of g(t). The time inte- 
gration only contributes a multiplicative constant and 
the gravitational wave spectrum is entirely determined 
by F(k) (the Fourier transform of the spatial structure 
of the source). From Eq. (fT0|) we find in this case the 
generic expression, y = (3{t — t*) 
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In Fig. [T] we show the second line of the above equation, 
namely (3 2 k 3 Ke[P s (k, k, k)], as a function of k/(3. This 
quantity determines the spectral shape of the GW spec- 
trum. We plot it for the three different forms of g(t) 
Eqs. (QUI, Qj]), (Ull), with \F(k)\ 2 given by Eq. (T^ and 
we choose two different velocities v = 1 and v = 0.01. 
Clearly, the shape of the GW spectrum is entirely de- 
termined by fc 3 |_F(fc)| 2 . The peak frequency corresponds 
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to k ~ R^ 1 = (3/v, the low frequency slope is k 3 and 
the high frequency one is 1/k. The different choices for 
g only slightly affect the amplitude but not the spectral 
shape which is entirely given by |-F(fc)| 2 . 
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FIG. 1: The function (3 2 k 3 Re[P s (k, k, k)] for the incoherent 
a function of k//3. The three curves correspond to g(t) 
given by (|16[1 . ((17} and (|18|) . and the velocities are v = 1.0 
(left curves) and v — 0.01 (right curves). 



B. Totally coherent sources 

Let us now consider the opposite extreme, when the 
source at different times is perfectly correlated, which we 
call totally coherent. We then have 

(n(k,t)n*(k,f)) - 

(2n) 3 S(k - k')y/P a (k,t,tWP a (k,t',t'), (21) 



so that 



P s (k,t,t') = \F(k)\ 2 \g(t)\\g(t')\ and 
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In this case the spectrum P s (k,k,k) is the product 
of the square of the space Fourier transform and the 
time Fourier transform of the source. Therefore, the k- 
dependence of the gravitational wave spectrum, namely 
the position of the peak and the power of decay at high 
frequency depends on the properties of the Fourier trans- 
form of g{t), denoted g(cj). Since the correlator has com- 
pact support in both, space and time, its Fourier trans- 
form is analytic in both k and u> hence typically starts 
with a constant. This plateau is expected to extend to 
the inverse of the duration of the source, (3, in frequency 
and to the inverse of the correlation scale, R^ 1 = (3/v in 
wavenumber. 

Since v < 1, the diagonal k = lu always leaves the 
plateau at k = oj = j3. Between (3 < k — lo < (3/v (the 



part of the diagonal between the horizontal and the verti- 
cal dashed lines in Fig. the function P s (k, k, k) decays 
with a power law depending on the assumptions on the 
continuity of g(t). For large lj, the Fourier transforms 
decay the faster the smoother the function is: we find 
the behavior uj~ 1 ,uj~ 2 and u>~ 3 for the three functions 
gi,g2 and g 3 , respectively defined in (fTrJ)) , (fTT|) and (fTSj) . 

For k > [3/v, P s {k,u),u)) decays even faster due to the 
additional suppression coming from the contribution of 
the spatial Fourier transform. In Fig. [2] we show schemat- 
ically the behavior of P s (k,uj,uj) in Fourier space. 

P (k, o, co) 



co = |3/v 



co = p 




k = p/v 



FIG. 2: The qualitative behavior of the function P s (k,ui,ui) is 
shown for the totally coherent case. The diagonal, P 3 (k, k, k) 
is also plotted. In the region uj < f3 and k < (3/v we expect 
a white noise spectrum of the anisotropic stress. For u) > (3 
and k > (3/v the spectrum is expected to decrease. Since 
the gravity wave spectrum only probes the diagonal ui = k, 
we expect, in the separable case with constant = (3/v a 
first change of slope at ui — k — (3 and a second at u) — k — 
(3/v. Whether the first or the second is the peak frequency 
depends on the space and time continuity and differentiability 
properties of P s {k,uj,u;). 

In Fig. [3l we plot the GW spectral shape 
(3 2 k 3 Re[P s (k, k, k)], for the coherent case, as a function 
of k/P, with F(k) from Eq. (T3J) and g(t) from Eqs. (IBJ 
to (|18[) . and for two choices of v. The plots confirm 
the qualitative expectations: for intermediate frequen- 
cies, (3 < k < (3/v, the slope of k 3 Ke[P s (k, k, k)] is linear 
in k if g(t) is discontinuous, it behaves like 1/k if g(t) 
is continuous, but the first derivative has discontinuities, 
and like 1/k 3 if g(t) is continuously differentiable once 
but the second derivative has discontinuities. It is inter- 
esting to note that only the behavior of the correlator 
close to the least differentiable points, i.e. the beginning 
and the end of the source is relevant for the behavior at 
large frequencies. 

For high frequencies k > (3/v we have the same be- 
havior discussed above, but multiplied by the decay of 
|F(/c)| 2 (which behaves as 1/k 4 ). These features are 
clearly seen in the plot with v — 0.01, where the interme- 
diate and high frequency regimes are well separated. It 
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is important to notice that the property of differentiabil- 
ity of g(t) influences the peak position, changing it from 
k = /3 to k = fl/v if the source is discontinuous in time. 




FIG. 3: The function /3 2 k 3 Re[P 3 (k, k, k)] for the coherent 
function of k//3: top panel, v = 1, bottom panel, 
v = 0.01. The three curves correspond to g(t) given by (|16[l . 
(|17|) and (|18|) . Notice the different peak positions for g\ with 
respect to the other two. 



C. Sources with a 'top hat' correlation function 

This case represents an intermediate possibility with 
respect to the two situations considered above: for a 
given wavenumber k, the source is correlated only if the 
time separation is sufficiently small. Given a parameter 
x c of order unity, the correlation is different from zero if 
\t — t'\ < x c /k. To realize this behavior we set 



thus 

P s (k,t,t') 



\F(k)\ 



P s {k,k,k) = |F(fc)| 2 Re 



g{t )\ 2 o(t'-t)e(^-(t>-t)) 

symmetric t <-> t ] 

t.+i 



dte^\g{t)Y 
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dt' e- tkt ' 



(24) 



In this case again, the GW spectrum bears no direct re- 
lation to the time Fourier transform of g(t), but it has a 
more involved behavior. In particular, if k is large, the 
upper bound of the second integral is always given by 
x c /k + t and we find 



P s (fc, k, fc) 



F(k)\ 



2 sin(:r c ) 



t, + j 



dt\g(t)\ 2 . (25) 



The remaining time integral only contributes as a con- 
stant. This is indeed what is shown in Fig. [4j where 
P 2 k 3 Re[P s (k, k, k)] is plotted as a function of k/f}, again 
with the same choices for F(k), g(t) and v as in the previ- 
ous examples. The situation is similar to the incoherent 
case, in particular the peak position is always k = [3 / 'v; a 
change in the slope from fc 3 to k 2 is observed when k > (3 
and approximation J25|) becomes relevant. 



10 



■H 10"' 

1 10"' 
^ 10"' 



10' 



""1 1 1 111 






- g!(t),V=1.0 






- g 2 (t), v=1.0 






- g 3 (t), v=1.0 






- g,(t), v=0.01 






- g 2 (t), v=0.01 






— g 3 (t), v=0.01 




if 





0.1 



10 

k/|3 



100 



1000 



FIG. 4: The function /3 2 k 3 Re[P s (k, fc, fc)] for the top hat case 
(with x c — 1), as a function of k//3. The three curves corre- 
spond to g(t) given by (|16|l . (I17|l and (| ISp and the velocities 
are v = 1.0 (left curves) and v — 0.01 (right curves). 



D. Stationary sources 



(n(M)n*(M')> 



(27T) 



-6(k-k')[P s (k,t,t) 



Q(t' - t)G - (*' - tfj + symmetric t <-> t' 



Although it seems contradictory to call a source "sta- 
tionary" which by definition is active only in a finite pe- 
riod of time, 1//3, this assumption has been considered 
in the literature [111 ], and so it is interesting to study it 
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also here. Furthermore, on timescales which are much 
shorter than the duration of the source, stationarity may 
be a viable approximation 2 . A stationary source is one 
where the unequal time correlator only depends on the 
time difference, 

(n(k,t)n*(k,i')) = (27r) 3 <5(k — k')P s (fe, t — t'), (26) 

and therefore 

P s (k,t,t') = \F(k)\ 2 \g(t-t>)\ 2 , (27) 

where the function g(t) now has a different meaning. In 
principle it is still a function with compact support, since 
the argument satisfies — 1/(3 < t — t' < 1 / 0, but continu- 
ity at the boundaries is no longer an issue. A straightfor- 
ward calculation gives for the spectrum (with t! — t = r) 

P s (k,k,k) = \F(k)\ 2 ( U + l) dt f' + ' 3 1 dre lkT \g( T )\ 2 . 
■It, Jt,~t 

(28) 

However, this is not the expression used in the literature, 
where instead, exploiting the stationarity, the above dou- 
ble integral is simplified to: 

P s (k,k,k) = ffi^ r dre lkT \g(r)\ 2 . (29) 

P J — oo 

This approximation holds if g is negligibly small for 
r < t* - t and t > U + 1/0 — t for all i* < t < t* + 1/0. 
Then the second integral can be extended to infinity 
and the first integral contributes just the duration of the 
source, 1/0. In the literature usually a Gaussian func- 
tion is chosen, g{r) = exp(— (r/3) 2 /2) (possibly with an 
extra k dependence, cf. [llj). In this case, the spectral 
function decays exponentially for k > while it behaves 
like |F(&;)| 2 for small k. 

To summarize this section, we found that the peak fre- 
quency is typically given by the correlation length of the 
source v/ (3 and therefore strongly depends on v, except 
under the coherent approximation when it is given by 
the characteristic time scale of the source 1/(3- Still, in 
this case, there is a kink at the frequency v/ (3. As for 
the high frequency part of the spectrum, it typically de- 
cays asl/fcocfc 3 |F(fc)| 2 , although in the coherent case, it 
strongly depends on the time structure of the anisotropic 
stress, especially its differentiability properties. 

IV. GW SPECTRUM FROM COLLIDING 
BUBBLES 

In this section we discuss the GW spectrum arising 
from bubble collisions. Especially, we want to compare 



2 Note that here we only consider stationarity of the source, and 
not of the induced gravitational waves as done in (ill . 



the results of Refs. Q and @ and comment on the dif- 
ferences. 

All the models in the last section have one feature in 
common: Due to the assumption of separability, Eq. (|14[) . 
the slope of the spectrum changes when the frequency 
surpasses the length scale of the problem, R^ 1 ~ (3/v. 
However, this feature is not seen in the GW spectra 
resulting from numerical simulations of bubble colli- 
sions |5j. Hence, one has to relax this assumption of 
separability to model the case of colliding bubbles cor- 
rectly. This was done in the analytic approach of Rcf. 6] 
where the length scale corresponds to the time-dependent 
bubble radius, R = v(t — t*). Nevertheless, this analytic 
approach leads to a distinct peak at k ~ R,^ 1 = /3/v, 
which is not seen in the numerical simulations. 

Of course, the case of colliding bubbles is quite special, 
and our goal is to find an analytic description that repro- 
duces most features found in the numerical simulations. 

We consider a bubble from a first order phase tran- 
sition which collides with a second bubble at time ti 
and equilibrates to a new, larger spherical bubble at 
time tf (or is absorbed by surrounding bubbles) and 
tf — U ^ ■ Let us assume that the tensor type 
anisotropic stress of this collision process is given by some 
function /„(x — x„, t — t n ) — /(y, r), where x n is the cen- 
ter of one of the bubbles which collide. We consider the 
function f n to be of compact support in both, space and 
time, continuous in time but with a kink at t — ti. This 
feature has been found inj| and it is confirmed in the re- 
cent simulations of Ref . [5| . The momentum density may 
be in the rapidly expanding bubble wall or it may also 
be in its interior. The tensor type (spin two) anisotropic 
stress is due to the fact that spherical symmetry is broken 
during the collision. 

Let us now constrain P s (k, fc, k) from this information. 
For simplicity, we suppress the tensor indices which are 
irrelevant for our considerations. The anisotropic stress 
power spectrum is given by 

(n(M)n*(k',t')>= (30) 

N N 

^^(e'^- k '^)/ n (k,t-g/;(k',i'-g) . 

n— 1 m—1 

This is the expression for the spatial Fourier trans- 
form of n(x, t) from N collision processes. Here f n is the 
Fourier transform of the tensor anisotropic stress from 
the n-th collision process which is centered at x„. We 
assume the center positions to be uncorrelated. There- 
fore 

(el ( k .x„-k'.x m)) = 2V- 1 S nm S(k - k') . (31) 

Here the volume V is included to take care of the di- 
mensions. This is a very reasonable assumption. First, 
the factor <5(k — k') is required by spatial homogene- 
ity. Besides, bubbles that are not in contact with each 
other should not lead to coherent effects and the cor- 
relation between overlapping bubbles is approximately 
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taken into account by the factor 2. This assumption of 
non-correlation then also insures that the total observed 
radiation reaches a constant value in the limit of large 
volumes. 

As we shall see now, only the density of bubbles, N/V 
enters in the physical result. We can write 



P s (k,t,t') 



P s {k,UJ,Uj') = 



P s (k,U!,0j) 



(2tt) 3 V 
2 

(2tt) 3 V 
2 

(2w) 3 V 

2 N 
J^fV 



: i(u)—u)')t ri 



(/ n (k, W )/*(k,w')> 



£<i/„(k,")i s 



l/(fc,")l 5 



(32) 



The function / in the last equation is the Fourier trans- 
form in space and time of a 'typical' bubble collision 
event. It is independent of the direction k because of 
the statistical average. Therefore, once we have deter- 
mined the Fourier transform of the anisotropic stress for 
a 'typical' collision process, we can just multiply it by 
N/V, the density of collision events, to obtain P s (k, k, k) 
and in turn the GW energy density spectrum. 

Notice that, because different bubbles are uncorrelated 
(c/. Eq. |3"T|) , the time Fourier transform enters in the 
spectral function (f3"2"|) . This feature is reproduced only 
in the coherent case, Eq. (f2"2"|) . Therefore, among the 
different cases discussed in the previous section and in 
Ref. @, only the coherent case can possibly reproduce 
the result obtained in Ref. |5|. 

Like in the models of the last section, the collision pro- 
cess has compact support in both, space and time, hence 
its Fourier transform is analytic in both k and lo and 
therefore typically starts with a plateau. This plateau 
is expected to extend to the inverse of the collision time 
scale, (3, in frequency and to the inverse of the typical 
bubble size, i?" 1 = j3/v in wavenumber, where v denotes 
the speed of the bubble wall. It can be deduced from the 
simple two bubble case Q , that (k, t) is continuous in 
time but its derivative is not: in particular, it has a kink 
at the initial time of action of the source. As a result, 
for frequencies k ~ /3, we expect the time Fourier trans- 
form to behave like in the coherent case of the last section 
Eq. (|22|) , in combination with the time-dependence given 
in Eq. (|TT)) (note that the term 'coherent' here refers only 
to the temporal behavior, while different bubbles are spa- 
tially uncorrelated). In particular, f{k, w) decays as 1/cj 2 
for large frequencies. The k-dependence of / is constant 
for k < i? -1 and is expected to decay for k > Hence 
the GW spectrum, which is proportional to k 3 \f(k, k)\ 2 , 
scales as k 3 for small frequencies, and beyond k ~ (3 it 
scales as at least up to k ~ R . Beyond we 
would expect it to decay faster than this behav- 

ior depends on the spatial dependence of the anisotropic 
stress. 



In the following we present a simple model that modi- 
fies the analytical model of Ref. [6| , in order to reproduce 
most of the features found in simulations of bubble col- 
lisions in the envelope approximation carried out in []|. 
These features are 

• For small wall velocities, the amplitude of the GW 
spectrum scales as v 3 and has a peak at a frequency 
k ~ (3. The peak position does not (or only very 
weakly) depend on the wall velocity [3j. 

• For large frequencies, the spectrum scales as k~ x , 
independent of the wall velocity even beyond 
k ~ R- 1 = 0/v. 

• For large wall velocities, the amplitude and peak 
frequency are slightly reduced [5[ (meaning the am- 
plitude grows slightly slower than v 3 ). 

We model the time dependence of the collision process 
by the function g2 (t) given in Eq. (|17|) which has the 
differentiability property we are looking for. The spa- 
tial Fourier transform might now be approximated by 
the expression (|15|) with R = v(t — t*). To recover the 
modeling used in Ref. @ we slightly modify the spatial 
Fourier transform to 



\f(k,t)\ 2 = R 3 



l + ( 



kR\2 
3 I 



i + (¥) 2 + (¥) 6 



(33) 



This is the result of Ref. [6[ for the spatial Fourier trans- 
form of the anisotropic stress, P s (k, t, t). Using the above 
expression for the bubble radius, it is easily seen that (|3"3"|) 
leads to a kink roughly at k ~ R^ 1 = f3/v, which is the 
size of the largest bubbles at the end of the transition. 
Because of the discontinuity of the anisotropic stress at 
the end of the transition, Ref. @ actually found a peak 
at R~ . However, the simulations only show a peak at /3 
but no peak nor a kink at f3 /v. The discrepancy between 
the two approaches is due to two different time evolutions 
for the correlation length. 

The simulations evaluate the time evolution of the por- 
tion of un-collided bubble wall. It is hereby assumed that 
the anisotropic stress is localized in a thin shell close to 
the bubble wall and that after the collision of neighbor- 
ing bubbles the stress vanishes inside the bubbles. Hence 
completely collided bubbles (whose walls are completely 
within neighboring bubbles) do not contribute to the 
anisotropic stress. Close to the end of the phase tran- 
sition the relevant length scale is then given by the di- 
mensions of the still un-collided bubble wall regions. The 
correlation length of the analytic model R(t) should then 
be replaced with this characteristic size in order to ap- 
proximate the simulation result. Indeed, once the transi- 
tion comes close to completion, even though the bubble 
sizes do grow, the typical size of colliding regions is ac- 
tually decreasing and tending to zero at the end of the 
phase transition (see Fig. [5]). Therefore, we replace R(t) 
by the size of a typical colliding region, which vanishes 
not only at the beginning but also at the end of the phase 
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transition. This reflects the fact that the source reaches 
a peak and eventually switches off. We model this by 
introducing a new characteristic length 



L(t) = -g 2 {t). 



(34) 



which we insert into ([33]) in the place of R and into the 
formula for the spectrum in the coherent approximation, 
Eq. (J2D) . However, when doing so we multiply by a factor 
L(t) 3 / 2 , and loose the property that f(k,t) should be 
C° but not C 1 at the endpoints of the transition (which 
would give us the correct slope). 

This problem can be fixed by arguing that the pre- 
factor R 3 , instead of being connected to the correlation 
length as in Q, actually just represents a volume factor 
connected to the un-collided bubble portion. Therefore, 
R 3 should rather be replaced by R 3 — ► L 2 AL, where 
AL ~ i?*e is a typical shell thickness. In this case the 
pre- factor becomes L(t)^/R^e. Here e < 1 is an arbitrary 
constant which is small in the thin wall approximation, 
the case considered in simulations. This (somewhat ar- 
bitrary) construction leads to 



f(k,t)=L(t)(j 



1/2 



1/2 



i + (^) 2 + (^) 6 



(35) 



where L = L(t) is given in Eq. 

The form of the resulting GW spectrum is shown in 
Fig. [6l This simple model leads to a velocity indepen- 
dent peak frequency and nicely reproduces the k~ 1 de- 
cay for large frequencies in accordance with the results 
from simulations 5]. In particular, there is no additional 
suppression at very large frequencies, k 3> f3/v. The os- 
cillatory behavior should vanish if averaged over several 
bubbles with slightly different nucleation times and sizes. 

Even though this model reproduces all qualitative fea- 
tures found in the simulations of bubble collisions in the 
envelope approximation, this analysis should not be un- 
derstood as a derivation, since some features of the spec- 
trum result from the judicious choice made in Eq. (f3"5)) . 
as we have argued above. For example, if we choose 
to replace L(t)y/R*e in the pre-factor by L(t) 3 / 2 , which 




±m t 



FIG. 5: Typical time evolution of the correlation length (cor- 
responding to the characteristic scale of the colliding region) 
and therefore of the source (anisotropic stress) generating 
the gravitational waves calculated in numerical simulations 
of bubble collisions [5| . 




FIG. 6: Qualitative behavior of the GW spectrum for the 
model of Eq. (|35[) that reproduces the results from numerical 
simulations of bubble collisions. 



seems more consistent, we obtain a.k 3 behavior for large 
k since the time-dependence is now C 1 , see Fig. [7] On the 




FIG. 7: The qualitative behavior of the GW spectrum is 
shown for the model in Eq. (|35[) . however replacing L(t)\/R,e 
byL(t) 3/2 . 

other hand, if we argue that the correlation scale should 
be the size of the largest bubbles R(t) — v(t — t*) as in 
, and we just fix the discontinuity problem by 



f(k,t) = g 2 (t)R 3/2 (t) 



1 + (M)3 



1/2 



i + (¥) 2 + (¥) 6 



(36) 



the spectrum has a peak at j3, a behavior between {3 
and v/f3, and an additional kink at v/j3 beyond which it 
decays more rapidly, see Fig. [5J 

Let us compare these findings with the results of [f|. 
First, we want to stress that bubble collisions in the ap- 
proach of Ref. [f| are modeled using Wick's theorem. The 
anisotropic stress correlator comes from the product of 
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FIG. 8: The qualitative behavior of the GW spectrum is 
shown for the model of Eq. 



four velocities, 



(n u (x)n im ( y )) 



(ui(x)Vj(x)t!|(y)i; m (y)) 
6m(x-y)&/(x-y) + 
£«( x - y)6"m(x - y) 



(37) 



where £i m (x — y) = («j(x)u m (y)) and the velocity cor- 
relator is non-vanishing only if x and y are in the same 
bubble. Thus, the above products are non-zero if x and 
y are in the same bubble with center, say z for the first 
factor and with center z' for the second factor. The 
probability that z = z' is vanishingly small. Hence x 
and y belong to two different bubbles which therefore 
must overlap. In this treatment, the scale associated with 
P s (k,t,t') is the size of the overlapping region which is 
of the order of the size of a typical bubble, R(t). The 
only difference between the approach followed in Ref. [6] 
and the Ansatz in Eq. ((36)) is the pre- factor gi(t). This 
pre-factor is nevertheless quite important since it ren- 
ders the anisotropic stress continuous in time. In Ref. [6] 
the anisotropic stress is an increasing function up to the 
end of the phase transition and is then abruptly set to 
zero. For the coherent case, this modifies the spectrum 
especially between (3 and (3/v, since the high frequency 
behavior of a discontinuous function goes like 1/fc. For 
this reason, in the totally coherent approximation which 
is relevant here, the gravitational wave power spectrum 
found in Ref. 6] grows like k 3 k~ 2 = k in the interval 
(3 < k < (3/v and decays only for k > (3/v, leading to a 
peak at i?" 1 = (3/v (cf. Fig. [3]). This behavior is seen 
only if the stresses build up dominantly towards the end 
of the transition, something that is not seen in the sim- 
ulations. Once this discontinuity is removed by e.g. mul- 
tiplying with (72 (t), the treatment proposed in Ref. @ 
turns into the spectrum shown in Fig [8j for which the 
peak is at (3 (cf. also Fig. [3]), but which still has a kink 
at R~ . This kink does not appear in the simulations. 

The simulations discussed in [B| show two major differ- 
ences with respect to the analytic modeling presented in 
Ref. @. First, in the simulations, the anisotropic stress 



correlator goes to zero at the end of the transition and 
does not decay abruptly as assumed in the analytical 
modeling. Furthermore, the typical scale which enters 
the spectrum is not the size of the bubbles as assumed in 
the analytical modeling, but it is the size of the not yet 
collided region of overlapping bubbles. This size starts 
small, reflecting the initial smallness of the bubbles, and 
goes to zero towards the end of the transition when most 
bubbles are nearly fully collided. 

In order to account for these differences, and recover 
the same result for the GW spectrum, we can modify 
the analytical model in [|| as explained above, leading to 
Eq. (|3"5|) . It is important to remark that the most obvi- 
ous treatment of the volume factor as L 3 (t) does not lead 
to the spectrum obtained in the simulations. To recover 
the results of the simulations, the shell thickness i?»e has 
to be introduced. This is not surprising, since the simu- 
lations are performed in the envelope approximation. 

Note also that the 1/fc behavior found in the simula- 
tions is very sensitive to the time differentiability of the 
spatial Fourier transform of an average collision event, 
the function f(k,t). This function has to be in C° but 
not in C 1 . The derivative has to have a jump (but not a 
divergence) at either t* or t* + l//3, and the non- vanishing 
slope of the right-side derivative at t* (or left-side at 
t„ + 1/(3) may not depend on k. This behavior can be 
modeled with Eq. (|35p where it is important not only 
that L(t) goes to zero in a continuous but non differen- 
tiable way at both ends of the phase transition, but also 
that the fc-dependence of f(k,t) vanishes at these times. 

Finally, let us estimate the constant of proportional- 
ity between the GW spectrum and the anisotropic stress 
|/(k, lo)\ 2 . The total, dimensionless anisotropic stress 
density for a typical bubble collision event is of the order 



/(x,t) 



Pv 



Prad ~\~ Pvac 



(38) 



where < k < 1 denotes the fraction of the latent 
heat that is transformed into kinetic bulk motion of the 
plasma and finally into anisotropic stress 0. For in- 
finitely thin bubbles in vacuum, where only the Higgs 
field plays a role, k = 1 and p rac i = 0. If the phase tran- 
sition happens in a thermal bath, Kp vac /p t ot — vp where 
Vf denotes the typical velocity of the thermal bath parti- 
cles resulting from the interaction with the bubble wall. 
The volume of a typical bubble is given by 



p3 V 

R* = Tio 'X T7- 

(3 3 N 



(39) 



When calculating the space Fourier transform of (|38|) we 
obtain a volume factor of the order of the size of the 
bubble. Furthermore, the time integration gives roughly 
a factor 1/(3 so that we have 



|/(k = 0,w = ' 



\ Ptot / (3 : 



(40) 



10 



Inserting this in the anisotropic stress power spectrum 
Eq. §2$ yields 



2 3 
vac \ ^ 



P S (0,0,0)~ K 2 (^ , 



(41) 



For the gravitational wave energy density given in 
Eq. PH|) this yields, together with the typical behavior in 
the wavenumber obtained from the numerical simulations 
and from the modeling leading to Fig [SJ 



gtr 



40 



rad 



2 ( Pvac \ ( 7~L* \ 3 



3tt 2 

(fc//3) 3 fc < fc pcak 
(/3/fc) fc pcak < fc. 



(42) 



For small wall velocities »Cl, the wavenumber of the 
peak is roughly constant fc pca k — /3, while a slight depen- 
dence on the velocity is observed in the simulation result 
for big velocities. 

Note, however, that the velocity v relating the charac- 
teristic scale and the characteristic time corre- 
sponds to the speed of the bubble wall only if the phase 
transition proceeds through detonations. In this case, 
v is anyway larger than the relativistic speed of sound, 
v > l/\/3. In the deflagration case, on the other hand, 
the speed of the bubble wall is subsonic, and the bubble is 
preceded by a shock wave in the symmetric phase (while 
the broken phase fluid is at rest). It is the collision of 
these shock waves that eventually leads to the generation 
of gravitational waves. Therefore, also in the deflagration 
case, the velocity relating the characteristic length and 
time scales of the problem is supersonic (since it corre- 
sponds to the front of a shock wave) Q. In summary, 
values of v smaller than the relativistic speed of sound 
are not realistic, and have been considered here just for 
illustrative purposes. This means that the difference in 
the peak position between the analytical and the simu- 
lation result, although conceptually relevant, is probably 
negligible from the point of view of observations. How- 
ever, it is important to notice that since the simulations 
are carried out in the envelope approximation, they can 
only model the detonation case and are valid under the 
assumption of supersonic velocities of the bubble wall. 



V. CONCLUSIONS 

In this paper we have discussed some general consid- 
erations which determine the spectrum of gravitational 
waves from a phase transition, or from some other short- 
lasting cosmological events which leads to the forma- 
tion of anisotropic stresses. A first, relatively known 
result is that the gravitational wave energy spectrum, 
dQ, gw /dlog(k) always grows like fc 3 on large scale, i.e. 
scales much larger than all scales in the problem. Fur- 
thermore we have seen that, if the unequal time correlator 
of the anisotropic stress is totally incoherent or coherent 



only over less than one wavelength, the time structure of 
the event does not affect the spectral shape, which is then 
entirely given by the spatial structure of the correlator. 
This situation changes if the source is close to totally co- 
herent. Then the spectrum changes at the characteristic 
time scale of the problem to turn from fc 3 to 

1. fc, if the anisotropic stress correlator is discontinu- 
ous (in time) at the beginning (or the end) of the 
source. 

2. fc -1 , if the anisotropic stress correlator is C° but 
the first derivative jumps at the beginning (or the 
end) of the source. 

3. fc~ 3 , if the anisotropic stress correlator is C 1 but 
the second derivative jumps at the beginning (or 
the end) of the source. 

These slope changes are realized if the jump height is 
independent of fc; otherwise, the result is more compli- 
cated. In case 1, an additional change of slope is needed 
at the typical spatial scale of the problem, for the total 
energy density to remain finite. Whether there is an ad- 
ditional change of slope in the other cases depends on 
the details. If the spatial structure does not have any 
intrinsic time dependence, i.e. in the separable case (I14p . 
this is certainly expected. However if the typical spatial 
scale of the problem depends on time this may affect the 
decay for large fc. 

Numerical simulations || indicate that for gravita- 
tional waves from colliding bubbles the second case above 
is realized, in such a way that there is no additional 
change of slope at higher frequencies. This is because the 
typical spatial scale of the source L(t) tends to zero also 
at the end of the phase transition. This result, which is 
at odds with the naive expectation that the typical scale 
would be the bubble size, R(t) (which tends to i?* = v / [3 
at the end of the transition) is quite important. It im- 
plies a mild 1/fc decay of the gravitational wave signal at 
high frequency which is most relevant for the detectabil- 
ity of the corresponding gravitational waves, e.g. from 
the electroweak phase transition (see Ref. Q). 

In the simulations, the typical spatial scale of the prob- 
lem is connected to the portions of un-collided bubble 
wall at a given time, and therefore goes to zero both 
at the beginning (when bubbles have not yet started to 
collide) and at the end of the phase transition. The sta- 
tistical average is then performed by averaging the GW 
emitted by a given realization over several directions. On 
the other hand, in the analytical approach followed in jfj, 
the characteristic randomness of the problem is assumed 
from the beginning. Therefore, what matters are corre- 
lation lengths, and the most obvious correlation length 
of the problem is given by the bubble size. This is so 
in the analytical approach of @ which models the 'over- 
lap' of bubbles, and not directly the 'collisions'. In con- 
trast, in the simulations, the bubble size does not appear 
as an important scale of the problem. Therefore, in or- 
der to recover the simulation results from the analytical 



11 



model, one needs to identify the size of portions of un- 
collided bubble wall as the relevant characteristic scale. 
Consequently, the 'volume' factor R 3 coming from the 
spatial Fourier transform (c/. Eq. ([55]) ) also has to be 
modified; however, this cannot be done simply by set- 
ting it to the un-collided bubble portion size cubed. To 
recover the simulation spectrum, the portions of uncol- 
lided bubble wall must only enter as a surface portion, 
while the thickness should be taken as an independent 
constant. This is to be understood in the context of the 
thin wall approximation used in the simulations. Taking 
into account a time dependent finite thickness of the shell 
of stress AL(t) would tame the kinks in Fig. [5] and also 
introduce an additional length scale in the GW spectrum. 
Consequently, using a finite wall thickness in the simula- 
tions would most probably lead to a steeper slope in the 
GW spectrum compared to the result in the thin wall ap- 
proximation for / > l/AL(ifi n ). Accounting for a finite 
wall thickness in the analytic approach [6j leads to the 



same result, but the steeper slope starts at / > 1/R(t& n ) 
(c/. Eq.Hil and Fig.®). 

It will be important to study the implications of these 
results for the production of gravitational waves from tur- 
bulence and from stochastic magnetic fields. In the first 
case, the typical spatial scale of the problem most prob- 
ably does not tend to zero at the end of the turbulent 
phase, while magnetic fields are likely to be long-lived 
and therefore have to be treated differently. 
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